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ABSTRACT 



The frequency response was determined for a second order servo- 
mechanism with backlash in the system. Frequency response for sinu- 
soidal inputs were obtained by digital computer solutions of the 
equations of motion of the system which had inertia and viscous 
friction on both sides of the gear train. System parameters included 
input amplitude, magnitude of backlash, damping ratio and resilience 
of the gear teeth. 

The Control Data Corporation 1604 Computer was utilized and 
Fortran programming was employed. 

Results show the effects of the amount of backlash and the 
gear material resilience on the output to input ratios. The exist- 
ence of sub-harmonic resonance was found with large values of back- 
lash . 

The authors express their gratitude to Dr. George J. Thaler 
of the Department of Electrical Engineering for his guidance and 
encouragement and to Professor Douglas G. Williams of the U. S, 

Naval Postgraduate School Computer Center for his generous assist- 
ance . 
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1. BACKGROUND AND INTRODUCTION 



Most servomechanisms have a mechanical connection between the 
driving components and the driven or positioning unit. Usually this 
connection is a gear train which, in most practical cases, has a 
finite amount of backlash in the gearing mechanism. This inherently 
renders the servo system non-linear and requires analysis techniques 
which are based on approximations o: which necessitate use of laborious 
graphical methods which may rely on simplifying assumptions. The 
backlash between the gear teeth cannot always be eliminated without 
possibly introducing other, perhaps less desirable effects such as 
coulomb friction and/or increased viscous friction. Manufacturing 
costs may also prohibit consideration of gears with no backlash. 

The relative distributions of viscous friction and inertia on 
each side of the gearing have direct effects on the behavior of the 
load. Graphical phase plane and analog simulation analysis techniques 
do not readily adapt to practical consideration and adjustment of 
these, distributions. Also, they cannot hanale the case of non- 
plastic impact of the gear teeth. This has led to the solution of 
this problem by the use of digital computers. 

Several investigations have been made of second order servo 
behavior using the digital computer at the U, S. Naval Postgraduate 
School. These include, limit cycle investigation by New for plastic 
impact (Ref. l) and by Anderson and Luckett for plastic and elastic 
impact (Ref. 2), and a transient response investigation for plastic 
impact by Andrews and Kelley (Ref. 3 ) . 
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This thesis extends the coverage of this problem to a steady- 
state frequency response invest igation, with both plastic and elastic 
impact, for a periodic (sinusoidal) input of varying frequency and 
amplitude. The general approach to the investigation was similar to 
that used by Anderson and Luckett, The equations of motion of the 
system were normalized so behavior would be independent of the un- 
damped natural frequency of the system. The parameters inertia, 
viscous friction, elastic ^bounce** coefficient and damping ratio 
were in the same ranges of values used by Anderson and Luckett so 
that a correlation of results could be made if desired. 

Computations were made on a CDC 1604 digital computer using a 
Fortran System program which had definite advantages over the machine 
language programs used in the previous investigations mentioned above. 
Solution of the system equations was based on the phase plane approach 
as to whether the system was operating combined or separated. Even 
though the system was normalized, the computer program was made gen- 
eral to permit extension of the values of the parameters if desirable. 

This thesis includes investigation of the effects of the amount 
of backlash, the amplitude of the periodic input and the amount of 
system damping for both plastic and elastic impact of the gear teeth. 
The only non-linearity considered was backlash. 

The total system was piecewise linear in that one linear differen- 
tial equation described the motion of the system when the driving motor 
was combined with the load and two different linear differential eq- 
uations described the separated motions of the motor and load. 
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Boundary conditions for the differential equations were based on phase 
plane slopes and on equations based on the laws of conservation of 
energy and momentum which were related to whether gear contact was 
plastic or elastic. 

The original intention was to determine the effects of different 
distributions of inertia and viscous friction between the driving unit 
(motor) and the driven unit (load). However, because of the great 
amount of digital computer operating time required, it was necessary 
to confine the investigation to one ratio of inertias and one ratio of 
frictions . 

Results from the digital computer calculations provided data from 
which the magnitude of the output and the phase relation between input 
and output could both be determined. This thesis is concerned only 
with the magnitude of the output as compared with amplitude of the in- 
put sinusoid. 

Special effects or phenomena such as jump resonance and sub- 
harmonic resonance, which are known to exist under certain conditions, 
were looked for in the course of the investigation. 

In the following sections the system equations are developed and 
the evolution of the computer program is described. A print of one 
version of the Fortran language computer program and the computer 
solution flow chart are included as appendices. 
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2. SYSTEM DESCRIPTION AND DEVELOPMENT OF EQUATIONS 



2 . 1 Description of System 

For the purpose of this investigation, a second order system was 



chosen. The block diagram of the system is as follows: 




Fig. 2-1 

This system is the same system that was reported by Anderson and 
Luckett in Ref. 2. Except for the non-linearity, the system is a 
simple second order servomechanism. There is only one non-linearity 
present, and that is backlash. The backlash is present in the gear 
train which separates the motor from the load. The effect of the 
backlash is to allow movement of the load and/or the motor inde- 
pendent of the movement of the other, so long as this movement does 
not exceed the amount of backlash. Also included as part of the gear 
train is the gear ratio; however, this has no special effect on the 
system other than is normally expected. 

Separated by the backlash is the motor inertia and friction and 
the load inertia and friction. The distribution of these quantities 
cannot be changed in the describing equations of the system due to 
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the backlash* In the linear case, the inertia and friction can be 
lumped on either side of the gear train without affecting the system 
operation, but this cannot be done when there is backlash present. 

The motion of the system can be described by a combination of 
three linear differential equations. When the gear teeth are touching, 
a single differential equation describes the motion. Once the gear 
teeth separate, one linear differential equation describes the motion 
of the motor, and another describes the motion of the load. Upon 
recombining of the gear teeth, the single differential equation again 
describes the motion. At the instant of separation or recombination, 
the appropriate differential equations have as initial conditions the 
last values of the equation which described the motion prior to the 
separation or recombination. The change over from one describing 
equation to the pair of equations and vice versa will be discussed 
further in the equation development section. 

2 .2 Assumptions 

Prior to developing the equations which describe the motion of the 
system under consideration, the assumptions which have been made will 
be presented. These assumptions are essentially those made in Ref. 2 
and are repeated for convenience but with some modification. 

The following assumptions were made: 

1. The gear teeth were initially in contact, and the initial con- 
ditions of the system were all equal to zero. 

2. Plastic deformation of the gear teeth during steady state con- 



5 



tact or on impact and any torsional deformations of driving 
shafts are negligible, 

3. The inertias of the gears and drive shafts are considered as 
part of the load or motor inertia depending on their attach- 
ment in the system. 

4. The law of conservation of energy is completely satisfied in 
the perfect elastic case, and the law of conservation of 
momentum is satisfied in the imperfect elastic cases includ- 
ing the plastic case. When only the law of conservation of 
momentum is required to be satisfied, the energy lost by the 
system is dissipated in the form of heat by infinitesimal 
deformations of the gear teeth. 

The gear teeth are in contact only instantaneously during im- 
pact for all cases except the plastic case, and the impulse 
torques of drive, friction, bearing supports, etc,, are zero 
during impact, 

6. The coefficient of restitution of the two opposing gear teeth 
is the same or can be described by an equivalent coefficient 
if the gear teeth are of unequal coefficients. 

7. Backlash is assumed to be equal at all points on the gear cir- 
cumference, and the backlash is measured at the output shaft. 

2 .3 Equation Development 

From Fig. 2-1, the error signal, E, is the driving signal. This 
is equal to: 
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E - o R - o c 



( 2 - 1 ) 



The torque developed by the motor, T^, is equal to: 

T M -KE -K (0r“0c) (2-2) 

Equations (2-1) and (2-2) are true when the load and motor are combined 
and when they are separated. 

The load torque, T^, is: 

T l - J L 0 C -h f L Q c (2_3) 

When the load and motor are combined, the load torque when referred to 
the motor side of the gear train is: 

T lm ~j°T l =/>(j L e t + f L 9c) (2-M 

The forward or positive direction of the system is defined to be 

when : 

0C =/>0 M (2 ' 5 > 

and the backward or negative direction when: 

0 C -y°0 M + A (2-6) 

By differentiating either equation (2-5) or (2-6) with respect to 
time yields the same results. These are: 
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0 C -f G m 



(2-7) 



and 

0 C V 5 0 M (2 - 8) 

Substituting equations, (2-7) and (2-8), into equation (2-4) 
yields after rearranging: 

T lm (j L + f L 0 M ) ■ (2 ' 9) 

The driving torque of the motor must equal the load torque referred to 
the motor side (combined system only) and the load torque of the motor 
itself . 

T M = J M 0 M + f M 0„ + T lm (2 - 10) 

Substitution of equations (2-1), (2-7)* (2-8), and (2-9) into (2-10) 
after simplifying gives: 



■A iVA) A , K /> n _ K/> n 

e ‘ + (j, " (Jm vU) 0r 



( 2 - 11 ) 



Rewriting equation (2-11) 

0 C + 2. J (jL> n 0 C 60 n 0, = 60 n 0 R (2-12) 

where 



8 



K/> 

Jl 



(2-13) 



6J 



z 



and 
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fjA +/J_L 

Jm V^Jl 



(2-14) 



Equation (2-12), therefore, defines the motion of the system when the 
gear teeth are combined. 

Whenever the gear teeth are separated, there is no driving force 
on the load, and its motion is defined by the following equation: 

0 L + -4- 0 L =0 (2-15) 

To be noted at this point is the change in notation for the load. 
In order to distinguish between combined and separated operation, the 
subscript C or L is used respectively for the output. 

Combined operation: 

Output position 
Output velocity 
Output acceleration - 

Separated operation: 

Output position 
Output velocity 
Output acceleration - 



©c 

0C 

©c 

e L 

©u 

e L 
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Although this convention in notation was adopted, it should be remembered 
that the same quantity is described by two different symbols, depending 
on whether the system is combined or not. 

By the gear teeth being separated, the load torque is not reflected 
back to the motor side of the gear train. Therefore, equation (2-10) 
reduces to: 

T = J A 4- f G (2-16) 

I M U N/l ‘ I NA S ' 



or by substituting equation (2-2) into ( 2 - 16 ) and after rearranging 
reduces to: 



8 



tA 



40 M tTe,^e 



j 



M 



J 






'C 



j, 



bA 



R 



(2-17) 



Since this investigation is directed towards the frequency response 
of the system, the forcing function is: 



— A sin cot ( 2 - 18 ) 

If all angular frequencies are divided by the undamped natural angular 
frequency, 0)^, this effectively is the same as setting 0)^ equal to 1, 

By this normalization, the angular frequency in equation (2-18) is 
really cd/od^ . Therefore, the results are applicable to any simple second 
order system having any natural frequency 0)^, if the forcing function 
angular frequency is multiplied by the o> n of that system. 

Equation (2-12) describes the motion of the system when the gear 
teeth are combined. Equations (2-15) and (2-17) describe the motion of 
the load and motor respectively when the gear teeth are separated. 
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Left to be determined now is the instant when separation or combination 



occurs • 

Use of phase plane analysis is necessary to determine the instant 
of separation when the gear teeth are combined* By using velocity and 
position coordinates for the phase plane, the slope of a trajectory is: 

de 

d 9 





and specifically for the combined system: 



M = d9 c _ _0£_ 

c da " a 



( 2 - 20 ) 



With proper manipulation of equation (2-12), the following obtains: 



N, = 



_ CUn 0R ~ COn 0 C ~ Z^SlOy\ 0 c 



a 



( 2 - 21 ) 



Using similar techniques to find the load alone slope on the phase 
trajectory from equation ( 2 - 15 ) yields: 









1 L 

J, 



(2-22) 



Assuming that the system is combined, it will stay combined as 
long as the magnitude of the load alone trajectory slope is greater than 
the magnitude of the combined trajectory slope. Symbolically, the 
system will remain combined so long as the following is true: 
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(2-23) 



N. 



> 



N 



c 



Separation occurs when the combined trajectory slope equals the 
load alone trajectory. This, therefore, is the criterion to be satis- 
fied for transferring the motion description of the system from 
equation (2-12) to equations ( 2 - 15 ) and ( 2 - 17 )* 

For the system to recombine after separation, the gear teeth must 
come into contact with such velocities that they do not bounce apart in 
the elastic case. To recombine in the plastic case, the gear teeth 
merely must come into contact. 

In the general case, two criterion must be satisfied for the load 
to recombine with the motor. These are: 

1. The gear teeth must touch. 

2. The velocities of the load and motor must be such that they 
are equal after impact. (Elastic case only) 

Unless equation (2-5) or (2-6) is true, the gear teeth cannot be 
in contact. Therefore, when either equation (2-5) or (2-6) is satis- 
fied, the gear teeth have come into contact. In making this test, 
however, it must be remembered that is actually used in place of 

0 C in equations (2-5) and (2-6). It is re-emphasized that these 
symbols actually represent the same quantity but are obtained by solving 
different equations . 

If the impact is plastic, the gear teeth will stay together. How- 
ever, in the elastic case, they do not necessarily have to remain in 
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contact. During impact, the law of conservation of momentum must be 



satisfied; and after impact, the velocity ot the load and motor must 
be the same for the teeth to remain in contact. By using the lav; of 
conservation of momentum, the respective velocities can be determined. 
In Ref. 2, the coefficient of restitution was defined. This is; 



e 



-Q l ) 

(yO0^ - 0 L j 



(2-24) 



or rearranging: 



_ 



&l 6 (p Qm ej 



(2- 24a) 



where the primed terms represent the respective velocities after im- 
pact. To satisfy the conservation of mementum law: 



j 



p 



- 0„J - (0 L - 0 L 



(2-25) 



must be true. To be noted, equation (2-25) is independent of the 
coefficient of restitution. 

• / 

Solving equation (2-24) for 0^, gives equation (2-24 a); and sub- 
stituting the result into equation (2-25), the following expression for 



e L 



obtains : 



e L = 



J 



lA 



(Jmv 2 j l 



/ > 0 M (e + i) + Q c ~ e 



(2-26) 
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If the gear teeth have come into contact, equations (2-26) and 

# / • / 

(2-2^A) can be solved for 0 L and 0^ respectively. Unless 
•i • / 

= Qm > the gear teeth will not remain together. 

In summary, assuming that the system starts in the positive direc- 
tion from zero with the gear teeth in contact, the operation is as 
follows : 

1. Equation (2-12) defines the motion of the system until equation 
( 2 - 23 ) is no longer true, 

2. Once equation (2-23) is no longer true, equations (2-15) and 
( 2 - 17 ) describe the motion, taking as their initial values 
the terminal values of (2-12). In general, the direction of 

the system must be determined so that t;he proper value for 

» 

e« can be used to initialize equation (2-17)* If 0c is 
positive, equation (2-5) applies. If 0 C is negative, 
equation (2-6) applies. In transferring from equation (2-12) 

to (2-15), = 0 C . 



3. During separated operation, recontact is determined from equa- 
tion (2-5 ) or (2-6) . 

h. In the event of recontact, equations (2-26) and (2-2^A) are 
used to determine if bounce occurs, (Required only in elas- 
tic case) 

5. Finally, when there is no bounce, description of the system 1 s 

motion reverts to equation (2-ll). This time, the terminal 

value of 0 L is used to initialize equation (2-12) along 
• / 

with 0c . 

6. The steps of operation continue this pattern. 
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2 .4 Parameters 



By referring to equations (2-12), (2-15), (2-17), and (2-24), the 
following parameters are available for varying: 

1. y 

2. iO n 

3. 0 R 
^ Jm 

5 • f m 

6 . J i_ 

T. f L 

8. e 

By setting 6J n equal to one, the system was made independent of 
its natural frequency. This normalization was discussed previously. 
Not shown in the preceeding list as a parameter is the gear ratio, 
yO , of the gear train. Since the gear ratio was considered to have 
no special effect on the system, its value was set equal to one. Also 
not shown is the magnitude of backlash. This was not used as a para- 
meter; however, it was assigned a value other than unity. 

The motor inertia and load inertia, and respectively, 
were combined such that with jO equal to one their sum was one. 

(i.e., J_, + yO 2 j = 1) This means that any combination of J w and 

7 M s L M 

can be used, so long as their sum is one. 

By solving equation (2-13) lor the system gain constant, K, the 
following is obtained. 




(2-27) 
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Since go is one, (j^ -f y° 2 J^) is one and yO is one; the gain con- 
stant, K, must also be one a This, therefore, was the value assigned to 
K in this investigation, 

Anothe-r variable was defined by combining the load and motor 

friction, f and f , similar tp the way the inertias were combined, 

JL M 

This new variable called total friction is defined as follows: 




(2-28) 



The friction ratio of the system is defined to be f /f . 

L JL 

Substituting the value of one for cjo ^, yO and (j^ + yO 2 J^) in 
equation (2-lU), the following results 

2 7 = f M + f L (2-29) 



or 



j 




(2-29A) 



Assigning a value to the damping factor, Cf , fixes the sum of ( f + f^). 

To determine the value of f or f T , a choice also must be made for the 

M L 

friction ratio. Choice of and the friction ratio, therefore, fixes 

f w and f Y . 

M L 

For this investigation, the imput, 0^ , was a sinusoid defined by 
equation (2-18). The angular frequency of oscillation, go, was effec- 
tively normalized by making co^ equal to one. However, the values 
assigned to go were varied over the desired range of the frequency 
response. There is also no restriction on the choice of the forcing 
function amplitude A. 



16 



Equation (2-24) defines the coefficient of restitution of the gear 
teeth. This definition imposes no restriction on its value; practi- 
cally however, it must be less than one and greater than zero. 

Having so many parameters available, it was necessary to restrict 
their variation. Listed below are the values assigned to the various 
parameters. It was felt that this choice would produce the most in- 
teresting results. The method used to vary parameters will be covered 
in Section 3* 
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TABLE II- I 



Values Used 
for 

Parameters 



PARAMETER 


SYMBOL 




VALUES USED 


Natural frequency 


00^ 




1.0 


Gear ratio 


A 




1.0 


Backlash 


A 




.1, .05 


Load inertia 


Ju 




.5 


Motor inertia 


vJfA 




.5 


System gain constant 


K 




1 .0 


Motor friction 


"Tm 


.12, 


.24, .36, .48, .60, .72, .96 


Load friction 


f L 


.08, 


.16, .24, .32, .40, .48, .64 


Friction ratio 


H 




.4 


Damping ratio 


$ 


.05, 


.1, .2, .3, .4, .5, .6, .8 


Forcing function frequency 


00 




.1 to io.o 


Forcing function amplitude 


A 




10, 1.0, .5, .095, .12 


Coefficient of restitution 


e 




0, .6, .8 
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3* DIGITAL COMPUTER PROGRAMMING AND METHOD OF SOLUTION 



3 . 1 Digital Computer Program 

All of the frequency response results were obtained by com- 
putations using the Control Data Corporation 1604 high-speed digital 
computer. The computer programming was done using the Fortran system 
with program input either by magnetic tape or by punched cards. Com- 
puter output was put onto magnetic tape and results were printed 
from the magnetic tape by use of an IBM 1401 high-speed line printer 
unit . 

The Fortran system was used because Fortran programming permitted 
inserting the equations of the servo system into the program in forms 
very close to their actual mathematical forms. The use of Fortran 
also allowed a great deal of flexibility in easily modifying the program, 
in changing system parameters, in inserting input data and in extract- 
ing the output results in a properly labeled decimal format. 

A print of a typical Fortran source program used in this inves- 
tigation is shown in Appendix B. This Fortran source program, includ- 
ing sub- routines, represents about 1600 machine language instructions. 

The system equations to be solved by the computer were developed 
in Section 2. The three differential equations describing the motion 
of the combined system and the motion of the separated motor and load 
were solved by subroutines which were parts of the overall Fortran 
program. These sub-routines were based on the Runge-Kutta numerical 
method for the solution of ordinary differential equations using a 
four stage process. Each of the Runge-Kutta sub- routines utilized its 
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own integral sub-program. 

Even though the # system being investigated was normal ized (g> =1) 

and a unity gear ratio was chosen ( yO = 1), the program was generalized 
to provide for changing 03^ and yO to any desired magnitudes. Similarly, 
it was possible to select any value of system gain K, however a value of 
K = 1 was used for all calculation for reasons discussed in Section 2. 

As shown in Section 2, the forcing function or input to the system 
in this frequency response investigation was = A sin cot. The 
amplitude A of the sinusoidal input could readily be set prior to any 
computer run. The sinusoidal input was obtained from the Fortran 
Reference Library sine function which is an integral part of the Fortran 
system. The input frequencies a> for each frequency response run were 
put into the program in an array of data. Each 0) was consecutively 
inserted into the calculations by the program until all desired ud s 
were completed. A knowledge of Fortran programming is necessary to 
thoroughly understand how this is accomplished. 

It was necessary to choose integration intervals which would 
yield accurate results for the Runge-Kutta solutions of the differen- 
tial equations yet which would not require more calculations than 
necessary. It was found that for values of ^f^O.l the following 
relationship was satisfactory for determining the intervals. 
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Integration interval 



0.001 

co 



( 3 - 1 ) 



A 



t = 



This resulted in 2000 x 2TT solutions per input cycle for all values 
of input frequency ao. The values of At were also put into the pro- 
gram in an array of data and each value was inserted into the calcu- 
lation along with its corresponding (JO. 

With combination of low 03 and Cf <0.1, it was necessary to re- 
duce At to a value less than obtained from equation (3-1) in order to 
maintain sufficient accuracy. 

Since the damping ratio y for any run was dependent on the rela- 
tive magnitudes of f^ and f^, these parameters (y , f^ and f^) were 
set prior to each separate run and were constant for that run. All 
runs were made with a constant ratio of 




f M +/ a 



Even though the distribution of inertia between motor and load 

was even for all calculations, the program provided for changing this 

distribution to 80 $> of the total inertia to either the motor or load 

side. This was built into the program and arrays of values for and 

J were put in as data similar to the arrays of 0) and At, 

L 

The amount of backlash between the gear teeth was fixed before 
eaqh computer run; in fact the value used for all but one run was 
A =0.1 radian. A value of A = 0.05 was used in a special run to 
check the effect of changing A and A proportionally. 
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Because of the large number of solutions per input cycle ( >OOCTT ) 
and since it was not necessary to record all of these results, the 
program was designed to print* the results for every one hundredth 
position computation. This provided sufficient results to obtain the 
frequency response of the system. 

Only the steady state response, after the system transients had 
effectively died out, was desired. Therefore a delay of five time con- 
stants (based on a similar linear system) was included in the program so 
outputs would not be printed* until this transient period had passed. 
The program mnenomic for this was DELAY. 

It was necessary to tell the program how long to compute before 
switching to the next value of go. This was done by causing the compu- 
tations to shift to the next o> when a predetermined amount of problem 
time had elapsed after the transient period. This time was based on a 
desired number of input cycles. The program mnemonic for this was 
designated CYC. 

3*2 Sequence of Solution of System Equations 

* _ , --t ' “ “* i 

The starting condition for each frequency response run was with 
the gear teeth in contact, i.e. the motor and load were combined, with 
this equation describing the system motion: 

* * * 2 _ 2. 

0 C + 2 y GOy, 0 C + 0 C rr 60 n A 5 I n cot ( 3 - 2 ) 

This differential equation was solved as each increment of sinusoidal 
input was added. The corresponding value of 0 C was calculated and the 
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combined phase plane slope was computed by 



N 



c 




(3-3) 



This slope was compared with the computed value of the load-alone 
slope, which was constant for the run; 



L Ju 

When N - N T the motor and the load were commencing to separate. If 
the system velocity 0 c was positive at separation the gear teeth were 
separating from the original contact position. 

if e c was negative, the teeth were separating from the opposite side of 
the backlash. 

After separation, the equation describing the motor motion was 



_ a 



(3-M 



0 M + + ~7~ 0 L - —j— A Sin CO i (3-5) 

'Jm Jm 

This differential equation was solved for each input increment. At 
the same problem time t, the following differential equation for the 
load alone was solved: 



0 l + 4“© l - 0 (3_6> 

From these solutions, the positions of the motor and load were 
compared in order to determine when the gear teeth had re- touched on 
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either side of the backlash, A condition to be satisfied for re- 



touching was 



°F 



0 L -/0m 



0 L 0 + A 



(3-7) 

(3-8) 



When the teeth did make contact, the load velocity after touching was 

* 

calculated from 



. / 



e 



L 



Jm 

J m"*’// 



/> S M 6 + e) + 





(3-9) 



and the motor was 



/ 





(3-io) 



The following equation was then solved to see if the teeth stayed 
combined 



Is “ 0 L — © ? (3-11) 

Where € was an arbitrarily small number which was used to determine 
when the velocities were almost equal in the case of elastic impact. 

When the teeth stayed combined, the program returned to the 
combined differential equation (3-2) and the complete sequence was 
repeated . 
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3 *3 Printing of Results . 



Solution of any of the three differential equations of the system 
produced the values of the respective position @ and velocity © . 
Since the position and velocity of the combined system were defined as 
the load or output position and velocity, these values of 0 c and@ c were 
printed as the outputs when the gear teeth were iu contact. When the 
motor and load were separated, the load alone was printed as the 
output, and the motor position and velocity were also printed for the 
same problem time t. 

For each print of system position and velocity, the problem time 
and the input position were also printed. Prints were made for each 
one hundredth solution of the differential equation or equations con- 
cerned. Special prints of problem time, position(s), velocity(s) and 
input were made just after the gear teeth separated or combined. 

A typical print of computer output is shown in Appendix C. This 
shows that output position and velocity (© r and@ c ) only were printed 
when the system was combined, while motor and load positions and 
velocities also were printed when the teeth were separated. 

Only five decimal places of results were printed because of 
limitations in the accuracy of plotting results* The computational 
accuracy was somewhat greater since the Runge-Rutta method was 
accurate to about ( At) ; , which was 10 ^ for the largest At used. 

A complete computer program flow chart and explanation is in- 
cluded in Appendix A* 
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3 ,U Computation Time 

Computer operating time to obtain a steady state output varied 
considerably for different values of and for the different input 
frequencies CJD . For example, to compute through a five time constant 
transient period, with a low damping ratio of = 0.1 (long time con- 
stant) and with an a) of 2.0 (small integration interval), required 
about 28 minutes of computation. On top of this, it took about one 
and three-quarter minutes to compute and print the outputs for each 
input cycle. From one to eight cycles of input (after the estimated 
transient period) were required to obtain a steady state peak output. 
For the other extreme of a high of 0.8 and a low 05 of 0.1, the 

total computation and print time required to obtain an output peak 
took only about two minutes. 

It is apparent that the investigation required a great amount 
of computer time. All computer runs, including test runs, took a 
total of about 83 hours of computer operation. 

An attempt was made to develop a computer program which would 
approximate the energy stored in the system for steady state at the 
first input frequency and would retain this ^stored energy 1 ' when cc 
was changed. It was hoped this would materially reduce the transient 
period in reaching steady state at the new 0). Test runs were promising 
but time did not permit refinement of this program to where it could 
be used in the investigation. 
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k. RESULTS AND DISCUSSION 



k . 1 Determining M from the Computer Output 

From standard frequency response nomenclature, the symbol M is 
defined to be the magnitude of the ratio of system output to system 
input. Symbolically this is: 



Appendix C shows a sample print out of data from the computer. 

One of the quantities printed is output position. To determine 0 C , 
it was necessary to pick off the value of output position at its maximum 
or minimum. In general, there were two peaks. One was a positive peak 
and the other a negative peak. These two peaks came about because the 
input was permitted to go through only one cycle once printing began. 

The magnitude of these peaks were not always identical. This indicated 
that the system had not reached steady state; however, in most cases, 
the difference between the peaks was less than plotting accuracy. The 
average of the magnitude of these peaks was taken to be 0 C . 

With values of the damping ratio greater than about .4, the method 



about 0.3 and lower, there began to appear an oscillation at the natural 
frequency of the system super- imposed on the output wave form. Its 




(*.i) 



where 




of obtaining 




as described above was used. At damping ratios of 
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presence could be expected when a limit cycle would be predicted by 
Ref. (2). 

Fig. 4- la shows this phenomenon. The damping ratio, ^ , was 
0.1 in this case. Shown in Fig. 4- lb is the output waveform of the 
same system with a higher damping ratio. ( = 0.5) It can be 

seen here that there is no super- imposed "limit cycle". To obtain 
0 C | for the system whose outpuf,. is shown in Fig. 4- la, it was nec- 
essary to take the mean of the super- imposed oscillation. In the case 
shown, this was not particularly difficult. However, as O) was in- 
creased towards the undamped natural frequency, this procedure became 
more difficult. Fig. 4-2a shows the output waveform for the same 
system at 0) = .6. From this figure, it can be seen that it is vir- 
tually impossible to establish the mean of the "limit cycle" oscilla- 
tions. Even if the output wave shape were available for a longer time 
period, it would still be virtually impossible to obtain an accurate 
mean value. To compare the output at the same forcing frequency, but 
with higher damping ratio. Fig. 4-2b was drawn. Here =0.5, and it 
can be seen that there is no difficulty picking the peak output value. 

Since the method of picking 1 0 C | varied, it was not possible to 
program the selection of | 0 C J . In some cases, a programmed method 
would have been satisfactory, but unless all the data were examined. 



the value of 



Qc 



thus obtained would be questionable, 
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4 .2 Frequency Response Variation with Damping Ratio 

Frequency response curves, for plastic impact of the gear teeth, 
were obtained for values of-^-of 10, 5 and 0.95> with damping ratio 
y varying from 0.1 to 0.8. These curves are shown in Figures 4-3, 

4-4 and 4-5 respectively. For comparison purposes, a similar family 
of curves for a comparable linear system is shown in Fig. 4-6. 

It was determined by test runs, where both A and A were varied 

proportionally, that the frequency response was identical as long as 

the ratio A was cdnstant. For example, the same frequency response 

resulted with A = 1.0 and A = 0.1 as with A = 2.0 and A = 0.2. This 

makes the frequency response results generally applicable to all 

simple second order servo systems for the same relative amounts of 

backlash since the system was normalized with cc =1. For this reason 

n 

the abscissae of the curves are labeled cjd/o) , and the curves can be 

n 7 

applied to any simple second order system by multiplying by 0)^ of the 
system . 

The resonant frequencies, for values of damping ratio that resulted 
in resonance ( ^ 0.6), were found to increase as the amount of back- 

lash was increased. This increase in resonant frequency was greatest 
for the higher values of y while it was almost insignificant at a 
y of 0.1. 

The peak values of the ratio of output position to input ampli- 
tude at resonance M^, increased above the linear peaks as the amount 
of backlash was increased. The increase was greatest at the lowest 
value of y . A plot showing this variation from the linear for dif- 
ferent amounts of backlash is shown in Fig. 4-7- 
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As might be expected, the effect of backlash is greatest at the 
resonance frequencies. The frequency response curves, for all amounts 
of backlash, are above the linear curve at values of 0) below the re- 
sonant frequencies. At frequencies greater than the undamped natural 
frequency and for values of !/) 0.4, the response with backlash drops 
off more rapidly than the linear and the degree of drop-off is greater 
as the backlash is increased. In the same frequency range, the response 
curve at low values of is slightly above the linear curve (ignoring 
sub-harmonic resonance which will be discussed later). This indicates 
that with heavy damping, backlash results in a more rapid attenuation 
of output at frequencies above the undamped natural frequency of the 
system, while with light damping, attenuation is slightly less than the 
linear up to 2 CD . 
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4.3 E ffect of Magnitude of Backlash on Frequency Response for Constant 

A more detailed examination of the effect of the amount of backlash 
in the system compared with the input amplitude was made for a constant 
damping ratio. The value of =0.5 was chosen since it represents a 
practical damping ratio which is less than optimum and results in a 
moderately oscillatory system. Also, no limit cycle exists at this 
value of y according to Anderson and Luckett, Ref. 2, It was pointed 

out earlier that as long as the ratio of is constant the frequency 

response is invariant. Therefore, is a logical system parameter. 

Figure 4-8 shows the frequency responses at = 0.5 for f\/[\ of 
10, 5 and 0.95 with plastic impact of the gear teeth. The linear re- 
sponse is also shown on the plot. It can be seen that the response 

deviated from the linear even when the magnitude of backlash was only 
one- tenth of the input amplitude. The deviation increased as the rela- 
tive magnitude of backlash was increased. The effect of backlash was 
to change the response curve to resemble the response for a lower damp- 
ing ratio. 

At frequencies above the undamped natural frequency, with 
A/a = 0.95, the response curve dropped below the linear curve and 
reached a minimum near o> = 2 . Up to this input frequency the ratio 
of output frequency to input frequency was one to one. With a slightly 
higher cjd, the ratio of output to input frequencies f Q /f ^ rapidly trans- 
sitioned to l/2 and the magnitude of the output abruptly increased. 

This was the first of a series of sub-harmonic resonance peaks which 
were encountered as 0) was increased further. The actual values of M 
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between sub-harmonic resonance peaks were not all definitive as will be 
shown in the discussion of the output wave shapes for this region of a). 

At some frequencies where input frequency was not an integer multiple 
of the output frequency , the output magnitude could not be defined 
because of the random, unsteady nature of the load position. This is 
shown by a dashed line in Fig, 4-8. 

4 .4 Sub-harmonic Resonance with Plastic Impact 

The first sub-harmonic resonance peak occurred at co = 2.2, (Fig. 4-8). 
At input frequencies of 2.0 and below, the outputs were steady with ampli- 
tudes decreasing with increasing cjd, and the ratios of output to input 
frequencies were l/l. The wave shapes of the output for a) = 1,8, 1.9 
and 2.0 are shown in Figs. 4-9 and 4-10. At cjd = 2 . 1 the output was 
very unsteady (Fig. 4- 10b) and a value for M could not be determined* 

At this input frequency, the ratio of output frequency to input frequency 
f^/f^ was not constant but fluctuated with an average of between l/l and 
1 / 2 . 

When cjd was increased to 2.2, the frequency ratio f / f^ became a 
steady l/2. The output wave shape became steady but unsymmetr ical, as 
shown in Fig. 4-lla. The peaks of this output wave occurred at the times 
of contact of the gear teeth. The reason for the unsymmetr ical wave 
shape was due to the difference in relative velocities of the motor and 
load at the time of contact. Referring to Fig. 4-lla, the higher 
(positive) peaks occurred as a result of the lower relative velocity 
at the instant prior to the preceding contact (smaller peaks). This 
smaller relative velocity imparted a higher velocity to the load which 
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caused a greater excursion.. The opposite was true in going from contact 
at the higher peaks (high relative velocity) to the smaller peaks* This 
phenomenon is believed to be related to the fact that there was plastic 
impact of the gear teeth. The value of M at the sub-harmonic resonant 
peak at a) = 2.2 in Fig. 4-8 came from the average of the high and low 
peaks . 

When CD was increased to 2.3 the output wave became erratic, the 

ratio f Q /f^ was transitioning from l/2 to l/ 3 and the exact value of M 

was indeterminant but was obviously less than at CD =2,2 (Fig. 4- lib). 

At CD = 2.4 the output wave (Fig. 4- 12a ) was similar to the wave 

shape at the previous sub-harmonic resonance peak of CD = 2.2. The wave 

was steady but unsymmetrical and the ratio of f^/f. was a fixed l/ 3 • 

The value of M was further reduced (averaged from peaks of output wave). 

Going to CD = 2.5 (Fig. 4- 12b) the wave shape was again unsteady 

but f /f. was a consistent l/3 and M began to increase. This increase 
o' 1 

in M continued for the next four CD* s . 

Wave shapes for CD* s 2.6, 2.7> 2.8 and 2.9 shown in Figs. 4-13 a nd 
4-14, show very convincingly a sub-harmonic resonance. At CD = 2.6, the 
general wave shape of the output during this resonance is first observed. 
For all four of these frequencies, the system output had the same wave 
shape. As the frequency was increased, the M value also increased. 

All four of these frequencies gave a f /f of l/3. Just as the output 
wave shape was very steady, the f^/f ^ was equally as steady at l/3« 

The sub-harmonic rise associated with these frequencies is shown on 
Fig. 4-8. Since the value of M for these frequencies was so well 
defined, this portion of the curve was not dashed as was required for 
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frequencies between 03-2.0 and 03 = 2.4. 

A further increase in frequency beyond 03-2*9 resulted in a drop- 
off in M. Fig. 4-15 shows the output wave shape for 03 = 3*0. This 
wave shape was taken for a period of over 50 time constants . As can be 
seen, there are portions of this wave which are quite steady and 
portions which are very unsteady. The steady portion had a f Q /f ^ of 
l/3* The M value taken from this steady portion was used to plot on 
Fig. 4-8. No M value could be obtained from the unsteady portion due 
to its inconsistent movement. The f /f during the unsteady portions 
was not defined either. It averaged between l/3 to l/4 . 

At frequencies greater than 3*0, the output motion was erratic. 

This can be seen in Fig. 4-l6 for 03 - 3*1 and 4,0 as well as in 

Fig. 4-17 for 03 » 5*0. Since M could not be deduced from these wave 

shapes, no points are plotted on Fig. 4-8 beyond 03 = 3*0* 

The f Q /f^ for 03 - 3*1> Fig, 4-l6a, was 1 / 3 . This may be mislead- 
ing because only a small amount of the output wave was obtained. At 
03 =s 4.0, Fig. 4- 16 b, the f / f was between l/4 and 1/5 . Fig. 4~17> 
for 03 = 5*0, shows the output wave shape for a considerable time. It 

can be seen here that the wave was not consistent. The f /f. for this 

o' 1 

03 varied widely. Depending on the portion of the curve examined, the 
f Q /f ^ was as small as l/6 and as large as l/4. 

Although no sub-harmonic resonance was found beyond 03 = 2,S>, it 
would appear that it does exist. In transitioning from one sub- 
harmonic frequency to another, the output wave was quite unsteady. 

When a definite output frequency was established, the M was also 
defined. Since the f Q /f ^ ratio changed from 1/3 at 03 = 2*9 to l/6 at 



U) = 5 a sub-harmonic resonance at l/4 and l/5 must have occurred between 
these two frequencies. 

From the limited data available, it appears that the system tends 
to resonate at odd sub-harmonics. Frequencies between CD = 2.0 and 
0) = 2.2 did not have a fixed f Q /f ^ , Then at CD » 2.2, a steady l/2 ratio 
was observed. 

During the build up to the third sub-harmonic resonance peak, a 
steady ratio of l/3 was present. This steady f Q /f^ over several values 
of frequency, leads one to believe the system would tend to seek an odd 
sub-harmonic frequency. Immediately past the sub-harmonic peaks, the 
system began seeking the next sub-harmonic frequency. The even sub- 
harmonics did not appear steady except at CD = 2.2. There was a con- 
siderable range, CD =• 2.4 to CD = 2.9> where f Q /f^ was a steady 1/3 add- 
ing support to the belief that the system would seek odd sub-harmonics. 
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4 .5 Effect of Elastic Impact 



Frequency response data were obtained, at a constant 0.^* for 

values of a/a equal to 10, 5 and 0*95 with e ~ 0.6 and e ~ 0.8. These 
parameters were the same as for the results discussed in Section 4.3* 
except for the values of the coefficient of restitution e. The response 
curves for e ~ 0.6 and e - 0.8 are shown in Figs. 4- 18 and 4-19. The 
values of M for all three values of a/a were found to be slightly 
less for both cases of elastic impact than for e =r 0 (Fig. 4-8) at 
frequencies near resonance. This indicates that elastic impact tends 
to reduce as compared with plastic impact. There was practically no 
difference between the response curves for e - 0.6 and e - 0*8 up to 
input frequencies of 03 = 2. 

With a/A = 0.95 sub-harmonic resonance began to occur at about 
03 = 2.8 for e = 0.6. Up to 03 ® 2.7 (Fig. 4-20a) the ratio f /f was 
a steady l/l and the output wave shapes were steady and symmetrical. 

At 03 = 2.8 the output waveshape became erratic (Fig. 4-21a) and f Q /f^ 
was transitioning from 11 to l/2. The response curve remained above 
the linear at all frequencies (Fig. 4-l8), unlike the equivalent curve 
for e = 0 (Fig. 4-8) which dropped below the linear between and the 
occurrence of sub-harmonic resonance. 

The response curve reached a minimum at 03 ~ 2.8. At the next 
test frequency 03 ~ 2.9* the output increased and the ratio f / f ^ was 
between l/2 and 1 / 3 . No even harmonic was found as there was for e = 0. 
The wave shape at 03 = 2.9 was also uneven with very inconsistent peaks 
(Fig. 4-20b ) . At 03 = 3.0 and 0 ) = 3.2, f / f ^ was still between l/2 and 
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1/3 the value of M was increasing with unsteady waveshapes (Figs. 4-21b 

and 4-22a). At CD = 3*4 (Fig. 4-22b) f q / f ^ became a consistent l/3 and 

M remained about the same, but the wave shape was still not steady. As 

0) was increased to 3*7 a nd to 4.0 (Fig. 4-23) the waves became symmetrical, 

f /f. remained l/3 and the value of M was decreasing. At the highest 
o 1 

test frequency CD = 4.3* M was still decreasing, the wave shape was be- 
coming unsymmetrical (Fig. 4-24) and f / f . was decreasing toward l/4. 

* 01 

In comparing the response curves of Fig. 4-8 (e ~ 0 and Fig. 4- 18 
(e = 0.6) for a/a = 0.95* it can be seen that the addition of bounce 
changed the curve and the sub-harmonic resonance pattern at frequencies 
above 0)^. Bounce eliminated the drop of the response curve below the 
linear response and apparently al o eliminated the even sub-harmonic 
resonance peak. The odd ( 1 / 3 ) sub-harmonic resonance peak was lower 
with bounce but was wider, i.e. the value of f Q /f^ =• 1/3 persisted over 
a wider range of frequencies. Sub-harmonic resonance started at a 
higher frequency when bounce was introduced. 

For e = 0.8 and a/a - 0*95* sub-harmonic resonance started to 
occur at about the same value of CD as for e = 0.6. The investigation 
did not go beyond this frequency for e = 0.8. 

In general, elastic impact of the gear teeth tended to reduce the 
effects of the nonlinearity for both the magnitude of the output and the 
amount of sub-harmonic resonance. 

4 .6 Response Characteristics with Very Low Damping Ratio 

In an effort to find jump resonance and the factors influencing 
its occurrence, some runs were made with a damping ratio, y , equal to 
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0,05. Considerable computer time was required for calculating the out- 
put at such a low damping ratio. Therefore, only a limited number of 
frequencies were investigated. 

The basic parameters for this investigation were: damping ratio 

equal to 0.05 > A/A equal to 1 „2, coefficient of restitution equal to 
zero. The frequencies used were between 0.5 and 1.8 with emphasis on 
the frequencies near the undamped natural frequency. Results of these 
runs did not show any jump resonance; however, some interesting things 
were shown. 

For a linear system, the resonant frequency is always less than the 
undamped natural frequency. (Resonant frequency is the frequency at 
which the resonant peak occurs on a frequency response.) Regardless 
of how small the damping ratio gets, the resonant frequency is always 
less than the undamped natural frequency; however, as y approaches 
zero the resonant frequency approaches but never exceeds the undamped 
natural frequency . 

With backlash present, this is not always so. Fig. U~25 shows a 
resonant peak at a frequency between 1 .01 and 1.02 which was greater 
than It has been shown previously that the effect of increasing 

backlash was similar to reducing the damping ratio. In addition to 
this effect, it appears that with a sufficiently large backlash and 
small enough y , the resonant frequency is increased more than the in- 
crease associated with a reduction in y alone. This may not be 
limited to low damping ratios. It may have occurred at higher damping 
ratios but was masked due to insufficient data being taken in the 
vicinity of the resonant peak. 
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The height of the resonant peak for this low is considerably 
higher than for the linear case. However, extrapolation of Fig* 4-7 
gives approximately the same value for this peak as shown on Fig. 4-25* 

At frequencies near the resonant frequency, between 0*95 a ^d 1*05, the 
output wave shapes were steady. The frequency ratio was also steady at 
one to one , Fig* 4-26 shows the output wave shape at a) = 0*95 and 
Fig* 4-27 shows it at a) - 1*0. As *can be seen, these are steady. 

The output below o> - 0*9 was rather difficult to interpret. By 
the previously discussed methods for obtaining M, it was not possible 
to determine this quantify accurately. Therefore, values of M for these 
frequencies were not plotted on Fig. 4-25 • 

Fig. 4-28 shows the wave shape for O) ~ 1.2. Even though the output 
wave shape is steady at this frequency, the input frequency is not an 
even multiple of the output frequency. This uneven relationship has 
been mentioned previously (Section 4.4) as an indication of the onset 
of sub-haxmonic resonance. Therefore, the low damping ratio decreases 
the frequency at which sub-harmonic resonance appears. Insufficient data 
are available to determine the exact location of the first sub-harmonic 
resonant peak, but it appears to occur at a frequency near ox - 1.6. The 
M at this peak is approximately 9> and the frequenc) ratio is even at 
one to two . 
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4 .7 Test of Program and Validity of Results 



For there to be complete confidence in the results of this investi- 
gation, some mention should be made of the steps taken to prove the com- 
puter program from which the results were obtained,, Also, some mention 
of the validity of the results is in order. It is the purpose of this 
section to discuss these subjects. 

Since there was no easy way that the system investigated could be 
simulated on an analog computer; and since little information is avail- 
able on such non- linear systems, it was necessary to test the program 
by approximation. First of all, the program was divided into two basic 
parts. These were, the combined system portion and the separated system 
portion. As integral parts of these two sections were the basic decisions 
which had to be made to transfer from one to the other. 

It was rather simple to test the combined section. This was done 
by setting the backlash, A , equal to zero. Doing this effectively 
eliminated the separated system portion and the decisions required to 
transfer into that portion. With the backlash equal to zero, the output 
should be the same as the output from a linear system. Comparison of 
the M obtained from the print out data, as described in section 4-1, 
with linear results obtained by conventional means were very favorable. 

The slight differences could be attributed to the computer solution 
still being in a transient state. The data from the computer solution 
was taken after five time constants had elapsed. The following table 
shows this comparison. 
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TABLE IV- 1 



O'- .6 A. 0 



03 


M, computer 


Mj linear 


.i 


.99714 


.99716 


.2 


.98810 


.98821 


.3 


.97172 


.97197 


.4 


.94693 


.94694 


• 5 


.91180 


.91192 



It was not possible to test the separated portion alone. There- 
fore, it was necessary to test this portion as part of the overall test 
of the program. For a linear system, the ratio of forcing function am- 
plitude to the amount of backlash, A,/A , is infinite since the backlash 
is zero. Therefore, if this infinite ratio is approximated by a non- 
linear system, the results should be very close to linear results. To 
approximate the linear system, an A,/a of 100 was chosen. The results 
of this test run were within 2.5 percent of the linear response at all 
frequencies less than the resonance frequency. At frequencies above 
resonance, the non-linear response was less than five percent different 
than the linear response out to the highest frequency used. The maximum 
percent difference occurred at the highest frequencies where small dif- 
ferences account for higher percentages. 

Since the program was capable of approximating a linear system 
while actually going through all of the calculations necessary for a non- 
linear system; it was reasonable to assume that with a greater amount 
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of backlash, smaller 



A, 
A ' 



the program would give valid results. 



Even though the time solution of the differential equations was 
considered accurate, the method of determining M may have introduced 
inaccuracies. For damping ratios greater than 0.3* the method of 
obtaining M was quite accurate. The only improvement would have been 
to calculate for a longer time, but this would have required more time 
than could be justified by the increased accuracy. 

At damping ratios less than 0*3j the interpretation of the output 
waveform, 0 C , was not so clear cut. At these damping ratios, the 

t » ft 

litait cycle was super- imposed on the output wave. At some frequencies, 

O) less than about 0.4, the method used to obtain M was quite good. 

(Mean of *limit cycle f oscillations ) However, for frequencies between 
about 0.4 and 0.7> the M was not so nicely defined. It was in this 
frequency range for low ifs that the results were not as accurate as 
they were otherwise. At all frequencies above an a) of about 0.7> the 
results wert quite good except as dicussed in the sub-harmonic section. 

Mention has been made that data was taken after five time constants 
based on a linear system. The wave shapes at earlier times were looked 
at to determine if this delay were required* From these, it was evident 
that the system was in a transient state. Even after five time constants, 
there was evidence of transients. However, the transient affect was 
quite small. Depending upon the parameters, the transient was virtually 
gone shortly after five time constants had elapsed. The only exception 
to this was discussed as part of sub-harmonic resonance. In the cases 
where steady state was reached, it was reached in less than six time 
constants „ 
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if ,8 Appli cation of Results and Extension of Method 



While there is limited direct application of the results of this 
investigation to servo design, there are some points that could prove 
helpful in the design of certain types of systems. If a system with 
backlash were to be operated with a periodic inputs Fig. would be 

useful in determining how much backlash could be tolerated for a given 
input amplitude to meet specified peak outputs « Also, the results 
showing the effect of gear resilience can be of benefit in choosing 
a gear material, since higher resilience tends to reduce the effect of 
the backlash o 

The general form of the computer program can be adapted to calcu- 
late any particular frequency response of similar second order systems 
with known system constants (K, J , J , A, etc*) for specified para- 

li 1 A 

meters. The program can be very easily modified to obtain transient 
response of second order systems with any desired sets of constants and 
parameters* Similarly, it can be readily applied to steady state (limit 
cycle ) investigations . 

As has been pointed out, the method of obtaining frequency response 
results in this investigation requires a great deal of digital computer 
time. If further investigations were to be based on the same general 
approach, especially at very low damping ratios, a technique of going 
from one frequency to another without having to go through a new tran- 
sient period each time, should be devised. Such a program was briefly 
attempted, as mentioned in Section 3* which transferred boundary con- 
ditions of the differential equations at the first frequency to the 



appropriate differential equations as initial conditions at the next 
frequency. This was an attempt to utilize the ^energy*’ stored in the 
system to more rapidly set up a steady state condition for the next 
frequency. It appears this concept may be the only way that would make 
feasible this digital method to more extensive investigations. It could 
result in showing jump resonance when it exists, where the program used 
herein apparently did not. 
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Fig. 4-13 
OUTPUT WAVE SHAPES 
5 AA = .95 e 




57 



Fig . 4-1 4 
OUTPUT WA\ E SHAPES 
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Fig. 4-22 
OUTPUT WAVE SHAPES 
5 = .5 A/a = .95 
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Fig. 4-23 
OUTPUT WAVE SHAPES 
S - *5 A/^ = .95 e = .6 
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5, CONCLUSIONS 



Based on the findings of this investigation, the following conclu- 
sions are made: 

1* The system frequency response remains unchanged if A and A 
are varied proportionally, with all other parameters constant * 

2* The presence of backlash in the system has an effect similar 
to lowering the damping ratio of the system* 

3* The effect of backlash is more pronounced at lower damping 
ratios* 

4* Increasing the backlash causes a more oscillatory response. 

5 * With a very low damping ratio and a small value of A/A , the 
first resonant peak can occur at a frequency greater than the undamped 
natural frequency . 

6* At sub-harmonic resonance the ratio of input frequency to output 
frequency is an integer* 

7* The output magnitude is erratic when the input frequency is 
not an integer multiple of the output frequency c 

8. The odd sub-harmonics are more pronounced; and with elastic 
impact, the even sub-harmonic does not appear* 

9* The lower the damping ratio, the lower the frequency at which 
sub-harmonic resonance will occur* 

10„ With all other parameters constant, sub-harmonic resonance 
will occur at a lower frequency with plastic impact than with elastic 
impact * 
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11. The effect of elastic impact is to reduce the effect of back- 
lash e 

12. With a coefficient of restitution of either 0.6 or 0.8, the 
frequency response of the system is essentially the same. 

13. Combinations of low damping ratio and high frequency require 
excessive computer time to obtain the frequency response. 

1U. When a limit cycle exists (Ref. 2), an oscillation at the 
undamped natural frequency of the system is superimposed on the output 



motion , 
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APPENDIX A - FLOW CHART OF COMPUTER OPERATION 
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EXPLANATION OF FLOW CHART NUMBERS 



1 . 



2 . 



3 . 



4 . 

5. 

6 . 

7 . 

8 . 



9> 

10 . 

11 . 

12 . 

13. 

14 . 
15 - 
16 . 



Set initial conditions and assign values to constants and parameters 
of system equations- i.e„, , f M , f L , A, A, e, K,yO , 

Problem time for delay of five time constants, number of input 
cycles desired after five time constants. 

*Print heading for results. 

*Pi int run number and values of e, Jf 3 and A. 

• * 

*Print column headings for Time, 0 , 0 , 0 , 0 and 0 . 

L I. ri M K 

Set first initial conditions for combined differential equation 1 
Time - 0, 0 = 0, 0 f 0 

Select values of O) and At, 

*Print values of a), At,, and J, * 

Solve differential equation for combined operation 



Q c + 2 y uo n 0 C 4- 0 C - cOn A si n cot 

Increment time by integration interval t — > t + At v 

Calculate value of input 0.,* 

K. 

Does time equal five time constants? 

Is backlash zero? 

Solve for , 2 z 

Q c = -2yo>„a c -u„Q c + u>„ A sin cut 

Calcu late 




Calculate 

n l = 

Does Nl c - N ? 
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17. 



18. 



19 . 



20 . 



21 . 



22 . 



23. 



24 . 

25. 

26. 



27. 



28. 



29. 



30. 



31. 



32. 



33. 



34. 

35. 

36. 



Count a time increment At, 

Gone through ICC At 1 s ? 

Set At count to zero. 

*Print values of Time, © , Q and 0 0 

Does time equal five time constants plus time for desired number 
of input cycles after transient period'? 

Set initial conditions for second and third differential equations* 
Does time equal five time constants? 

*Print values of Time, 0 , 0 and 0^ „ 

Are gear teeth separating from original contact position plus A? 



Calculate 

Calculate 





e L -A 

JL 

f 



Solve differential equation for motor alone 

0 M + T + ^ = ^ A sin wt 

vJm 

Solve differential equation for load alone 

X® L = 0 

Increment time t — * t +- At. 

Calculate value ot input 0^* 

Does time equal five time constants? 

Count a time increment At. 

Gone through 100 At' s? 

Set At count to zero. 

*Print values of Time, 0 t> 0,, 0.,* and 0^ . 

* 1/ L y lr M R 
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37. 



Does time equal five time constants plus time for desired number 



38. 



39. 

40. 



41. 



42 . 

43 . 

44 . 

45. 

46. 

47. 

48. 



of input cycles after transient period? 

Are gear teeth touching in initial contact position? 

Are gear teeth touching in initial position plus A? 

• / 

Calculate load bounce velocity @ L 

Jm 



©l 






e) + G l j^ L - 



Calculate motor bounce velocity 0 

a' _ 8 l - eP 6 m + 6 Ql 

yO 



Do gear teeth stay together? 

Set initial conditions for combined differential equation. 

Does time equal five time constants? 

• * 

*Print values of Time, 0 T , @_ . 0. © and 0 

L 1/ M M R 

Go to combined differential equation* 

Completed calculation for all values of cu? 

Stop 0 



* 



* Print* 1 refers to patting the indicated items or their values onto 
magnetic tape for later printing. 
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APPENDIX B 

COMPUTER FORTRAN SOURCE PROGRAM 



1 

23 

24 
50 



PROGRAM BLASH 

DIMENSION DM ( 3 ) , 01(3), W(23), TINT(28), Yl(2), Y2(2), Y3(2) 
COMMON ZETA, WN, A, OMEGA, FL, DLOAD, FM, DMOTOR, GC , THEl 
ZETA = 0.5 

Ife = 8 : 6 , 

DELAY = 10.0 
A = 1.0 
BL = 0.1 
E = 0.0 
CYC = 6.3 
GC = 1.0 
RHO = 1.0 
WN = 1.0 

EPS = • 45474 735E- 1 2 
FORMAT (3F10.I) 



READ INPUT TAPE 


7, 


1 ( DL ( J ) , 


J = 


1, 


3) 


READ INPUT TAPE 
FORMAT (7F10.2) 


7, 


1 ( DM ( J ) , 


J = 


1 f 


3) 


READ INPUT TAPE 
FORMAT (6F12.10) 


7, 


23 (W(I), 


I = 


1, 


28) 



READ INPUT TAPE 7, 24 (TINT(I), I = 1, 28) 

FORMAT ( 58H 1 RE SULTS OF NONLINEAR SERVO WITH BACKLASH GRAHAM / LL 



1 0YO ) 

WRITE OUTPUT TAPE 8, 50 
2 OFORMAT ( 23H0 PRODUC T I ON RUN 19 E = , 

1F3.1, 3X, 7HZETA = , F4.2, 3X, 4HA = , F6.3) 

WRITE OUTPUT TAPE 8, 2, E. ZETA, A 

3 OFORMAT (5X, 4HTIME, I2X, 9H0UT POSIT, 13X, 7H0UT VEL, 9X, 1IHM0T0R 
1 POSIT, 11X, 9HM0T0R VEL, 11X, 9HINPUT POS) 

WRITE OUTPUT TAPE 8, 3 
DO 19 J = 2,2 
DO 19 I = 8, 9 
PAUSE 

DELT = T I NT( I ) 

OMEGA = W( I ) 

DLOAD = DL ( J ) 

DMOTOR = 0M( J ) 

FINI = CYC/Wl I ) 

T = 0.0 
' Y 1 ( 1 ) = 0.0 
Y 1 ( 2 ) = 0.0 

22 OFORMAT ( 8HOOMEGA = ,F5.1, 5X, l 1 HINTERVAL = ,F12.10, 5X, 

1 15HL0AD INERTIA =• ,F2. 1) 

WRITE OUTPUT TAPE 8, 22, W(I), TINT(I), DL ( J 1 
M = 0 

5 CALL RKUTTA1 (2, T, Y1 , DELT) 

T = T + DELT 

FORCE = A*SINF(OMEGA*T) 

THEC = Yl( 1) 

THECD = Y1 (2 ) 

I F ( T - DELAY) 26, 26, 30 
30 M = M+ 1 

IF (M - 100) 26, 27, 27 
27 M = 0 

WRITE OUTPUT TAPE 8, 4, T, THEC, THECD, FORCE 

4 FORMAT (F10.3, 2F20.5, 40X, F20.5) 

IF (T - DELAY - FINI) 26, 19, 19 

26 IF (BL) 5, 5. 21 

21 THECDD = -2 • 0* ZETA* WN* THECD - WN*WN*THEC + WN*WN*A*S INF ( W ( I ) *T ) 
SLOPEC = THECDD/ THECD 
SLOPEL = -FL/DKJ) 

IF (SLOPEC - SLOPEL) 6, 6, 5 

6 THEL = THEC 
Y3( 1 ) * THEC 
Y2 ( 2 ) = THECD 
Y3(2) = THECO 

IF ( T - DELAY) 16, 16, 33 

33 WRITE OUTPUT TAPE 8, 4, T, THEC, THECD, FORCE 
16 IF ( Y3( 2) ) 7, 8, 8 

7 Y2( 1 ) « ( Y3( 1 ) - BL)/ RHO 
GO TO 9 
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8 Y2( l ) = Y3( I )/RHO 

9 CALL RKUTTA2 (2, T, Y2 , OELT ) 

THEM = Y2( 1) 

THEMD = Y2 ( 2 ) 

CALL RKUTTA3 (2, T, Y3, OELT) 

THEL = Y3( 1 ) 

THELD = Y 3 ( 2 ) 

T = T + OELT 

IP"?? - 0lH?5‘SSf e }s:'ji 

31 M = M + 1 

IF (M - 100) 28, 29, 29 
29 M = 0 

10 FORMAT (F10.3, 5F20. 5 ) 

WRITE OUTPUT TAPE 8, 10, T, THEL, THELD, THEM, THEMD, FORCE 

IF ( T - DELAY - FINI ) 28, IV, 19 
28 IF (THEL - RHO* I HEM ) 13, 13, 11 

11 IF (THEL - RHO* THEM - BL ) 12, 13, 13 

12 CO TO 9 

IF (E) 17. 17, 13 

13 - BOLVEL = ( DM ( J ) / ( DM ( J ) + RHO * RHO * 0 L ( J ) ) ) *(RHO * THEMD 

1 * ( 1-0 + E ) f THLLD * ((RHO * RHO * DL ( J ) ) /DM(J) - E ) ) 

HOMVEL = (BOLVEL - E*RHO * THEMD + E*THELD )/RHO 
THELD = BOLVEL 
THEMD = HOMVEL 
Y2 ( 2 ) = THEMD 
Y3 ( 2 ) = THELD 

DIFVEL = RHO*THEMD - THELD 
IF (DIFVEL) 14, 15, 15 

14 DIFVEL = -DIFVEL 

15 IF (DIFVEL - EPS) 17, 17, 16 

17 Y 1 ( 1 ) = THEL 
Y 1 ( 2 ) = THELD 

35 IF ( T - DELAY) 34, 34, 1 8 

18 WRITE OUTPUT TAPE 8, 10, T, THEL, THELD, THEM, THEMD, FORCE 

34 GO TO 5 

19 CONTINUE 
END FILE 8 
STOP 
END 

SUBROUTINE RKUT I A 1 ( NUMBE R , XVAR , YVARS , STEP ) 

DIMENSION YVARS ( 30),AK(4 ,30) , D Y ( 30 ) , YC ( 30 ) , C ( 4 ) 

COMMON ZETA, WN , A, OMEGA, FL, DLOAO , FM, OMOTOR, GC , THEL 
' C( 1 ) = 0.0 
C ( 2 ) = 0.5 
C( 3) = 0.5 
C ( 4 ) = 1.0 
DO 1 1=1,4 
XC=XVAR + C ( I )*STEP 
DO 2 J=1, NUMBER 

2 YC ( J ) =YVARS(J) + C( I )*AK( 1-1, J) 

CALL DERIVl ( XC, YC,DY) 

DO 1 J=1, NUMBER 
1 AK ( I , J ) = STEP * DY( J ) 

DO 3 J = 1, NUMBER 

3 YVARS(J) =YVARS ( J ) + ( AK ( 1 , J ) +2 . *AK ( 2 , J ) +2. * AK ( 3 , J ) + AK ( 4 , J ) ) /6 . 
RETURN 

END 

SUBROUTINE DERIV 1 (1, Y, DY) 

DIMENSION Y( 2 ) , DY(2), W ( 28 ) 

COMMON ZETA, WN, A, OMEGA, FL, DLOAO, FM, DMOTOR, GC, THEL 
DY ( 1 ) = Y( 2) 

DY ( 2 ) = WN*WN*A*SINF(OMEGA*T) - WN*WN*Y(1) - ( 2. 0*Z ETA* WN ) *Y ( 2 ) 

RETURN 

END 

SUBROUTINE RKUT TA2 ( NUMBE R , XVAR , YVARS, STEP ) 

DIMENSION YVARS ( 30) ,AK(4,30) ,DY ( 30 ) , YC ( 30 ) , C ( 4 ) 

COMMON ZETA, WN, A, OMEGA, FL, DLOAO, FM, DMOTOR, GC, THEL 
C( 1 )=0.0 
C( 2) =0.5 
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C ( 3 ) =0. 5 

C(4 ) = 1 .0 

DO 1 1=1,4 

XC = X VAR + C( I ) * S TEP 

DO 2 J = 1 , NUMBER 

YC ( J ) =YVARS ( J ) + C ( l ) *AK ( I- I , J ) 

CALL DERIV?( XC,YC»DY) 

1 BE R 



YVARS(J) =YVAR S ( J ) + ( AK ( 1 , J ) + 2. * AK ( 2 , J ) +2. * AK ( 3 , J ) * AK ( 4 , J ) ) /6 . 

RETURN 

END 



L ALL UfcK IV/I a 

00,1 ‘It'jtNyyp 



SUBROUTINE DERIV 2 (1, Y, DY) 

DIMENSION Y( 2 ) , DY(2), W(28), DM(3) 

COMMON ZETA, WN, A, OMEGA, EL, OLOAD, FK, DMOTOR , GC, THEL 
DY ( I ) = Y ( 2 ) 

ODY ( 2 ) = (GC / DMOTOR) * A* S INF (OMEGA* T ) - (FM / DM0T0R)*Y(2) - 
I ( GC / DMOTOR )*THEL 
RETURN 
ENO 

SUBROUTINE RKUTT A3 ( NUMBE R , X V AR , YVAR S , STE P ) 

DIMENSION YVARSt 30 ) » AK (4 , 30 ) ,DY(30) ,YC(30) ,C(4) 

COMMON ZETA, WN, A, OMEGA, FL, DLOAD , FM, DMOTOR, GC, THEL 

C( 1 )=0.0 

C( 2 ) =0. 5 

C ( 3 ) =0. 5 

C ( 4 ) = 1 . 0 

00 I 1=1,4 

XC=XVAR + C( I ) *STEP 

DO 2 J=l, NUMBER 

> YC ( J ) =Y VARS ( J ) + C( I )*AK( I- 1, J) 

CALL DERIV 3 ( XC , YC , DY ) 

DO 1 J=l, NUMBER 
I AK { I , J ) = STEP * DY(J) 

DO 3 J = 1, NUMBER 

i YVARS(J) = YVAR S ( J ) + ( AK ( 1 , J ) + 2. * AK ( 2 , J ) +2. * AK ( 3 , J ) + AK ( 4 , J ) ) /6 , 
RETURN 
END 



SUBROUT 
^ DIMENSI 
COMMON 
DY( 1) = 
DY ( 2 ) = 
RETURN 
ENO 
END 

.2 
.8 
. 1 
.8 

2.5 

9.0 

.01 

.0014285714 

.0005882353 

.0001666667 

.0000333333 



INE DERIV 3 (T, Y, DY) 

ON Y(2), OY ( 2 ) , W ( 28 ) , DL(3) 
ZETA, WN, A, OMEGA, EL, DLOAD, 
Y ( 2 ) 

— (FL /DLOAD ) * Y{ 2 ) 



FM, DMOTOR, GC , THEL 



.5 

.5 

.2 

.9 

3.0 

10.0 

.005 

.00125 

.0005 

.0001428571 

.000025 



.8 
.2 
.3 

1.0 

4.0 
20.0 
.0033333333 
.001 1 1 1 1 1 1 1 
.0004 
.0001 25 
.00002 



.4 .5 

1.2 1.4 

5.0 6.0 

30.0 40.0 

.0025 

.001 

.0003333333 
.000 1111111 
.0000166667 



.6 

1.7 

7.0 

50.0 

.002 

.0008333333 

.00025 

.0001 



.7 

2.0 

8.0 

60.0 

.00 16666667 
.0007 142857 
.0002 
.00005 
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APPENDIX C 



TYPICAL PRINT 
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